| Section | Description | Status |
|---|---|---|
| Data summary | Summary of data | |
| Noise | General accounting/ cleaning the data to exclude noise | |
| Coverage | Plots of visual field coverage | |
| R2 | Variance explained plots | |
| Sigma | Sigma plots | |
| Eccentricity | Eccentricity plots | |
| Eccentricity x Size | Relationship between size and eccentricity | |
| Polar angle | Polar angle plots | |
| Noise | Plot that shows fit of the noise model | |
| Percent Signal Change | Gif of psc change over time |
library(ggplot2)
library(ggforce)
library(mixtools)
## mixtools package, version 1.2.0, Released 2020-02-05
## This package is based upon work supported by the National Science Foundation under Grant No. SES-0518772.
library(cowplot)
##
## ********************************************************
## Note: As of version 1.0.0, cowplot does not change the
## default ggplot2 theme anymore. To recover the previous
## behavior, execute:
## theme_set(theme_cowplot())
## ********************************************************
library(kableExtra)
library(stringr)
library(pracma)
library(viridis)
## Loading required package: viridisLite
Below is the dataset:
| X.1 | R2 | sigma | isigma | ecc | iecc | ang | iang | X | Y | ROI | idx | Hemisphere |
|---|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 27.695471 | 10.0000000 | 0 | 4.8539120 | 0 | 303.8647 | 0 | -4.0304727 | 2.7047647 | 0 | 0 | 1 |
| 1 | 56.418604 | 2.4727238 | 0 | 0.6329298 | 0 | 163.5531 | 0 | 0.1791997 | -0.6070318 | 0 | 1 | 1 |
| 2 | 4.057080 | 0.5649345 | 0 | 30.4843764 | 0 | 271.5453 | 0 | -30.4732891 | 0.8221032 | 0 | 2 | 1 |
| 3 | 1.943642 | 1.1748822 | 0 | 50.3402546 | 0 | 173.3346 | 0 | 5.8430497 | -50.0000000 | 0 | 3 | 1 |
| 4 | 44.977633 | 10.0000000 | 0 | 28.1751645 | 0 | 326.0413 | 0 | -15.7385219 | 23.3696132 | 0 | 4 | 1 |
| 5 | 1.031232 | 1.5382633 | 0 | 52.3031763 | 0 | 107.0665 | 0 | 49.9999999 | -15.3499922 | 0 | 5 | 1 |
Also define the images, and a function for making viridis colorbars (this code block should be hidden).
Set up some constants.
nROIS=length(unique(myfile$ROI))
facetval=floor(sqrt(nROIS))
recipe=list(scale_fill_manual(values=c("steelblue2","springgreen3")),theme_classic(),theme(axis.text.x = element_text(angle = 90)))
ROILABS=c("None","V1","V2","V3","hV4","VO1","VO2","LO1","LO2","TO1","TO2","V3b","V3a")
myfile$ROI=factor(myfile$ROI,levels=c(0:12),labels=ROILABS)
To try and parse between signal and noise, I fit a gaussian mixture model to the distribution of R2 values.
library(mixtools)
out<-normalmixEM(myfile$R2[!is.na(myfile$R2)],k=2)
## number of iterations= 49
Use this to define a threshold, whereby there is a low probability that the voxel belongs to the noise pool. Plot this.
thresh=out$mu[1]+(1.96*out$sigma[1])
myfile$noisesep=factor(ifelse(myfile$R2<thresh,0,1),levels=c(0,1),labels=c("noise","signal"))
varex=ggplot(data = myfile,aes(x = R2))+geom_histogram(color="black",aes(fill=noisesep),binwidth=.005, position="identity")+geom_vline(xintercept=thresh)+recipe
varex
## Warning: Removed 66365 rows containing non-finite values (stat_bin).
Now view this as a function of ROI
varexROI=ggplot(data = myfile[myfile$ROI!="None",],aes(x = R2))+geom_histogram(color="black",aes(fill=noisesep),binwidth=.01, position="identity")+geom_vline(xintercept=thresh)+facet_wrap(~ROI,nrow=facetval)+recipe
varexROI
Now, additionally determine the voxels that are tuned to the screen.
myfile$onscreen=rep(0,nrow(myfile))
myfile$onscreen=ifelse(myfile$X>-10 & myfile$X<10 & myfile$Y>-10 & myfile$Y<10,1,0)
Show position of prfs in relation to the screen.
ggplot(data=myfile[myfile$ROI!="None",],aes(x=X,y=Y))+geom_rect(xmin=-10,xmax=10,ymin=-10,ymax=10)+geom_point(aes(colour=onscreen))+recipe
ggplot(data=myfile[myfile$ROI!="None",],aes(x=X,y=Y))+geom_rect(xmin=-10,xmax=10,ymin=-10,ymax=10)+geom_point(aes(colour=ang))+recipe+scale_colour_viridis(option='plasma')
ggplot(data=myfile[myfile$ROI!="None",],aes(x=X,y=Y))+geom_rect(xmin=-10,xmax=10,ymin=-10,ymax=10)+geom_point(aes(colour=noisesep))+recipe
Before anything else, do some general accounting to determine the number of noisy/offscreen voxels in each ROI.
| N vertices | PropNoise | PropOffscreen | Mrsquared | |
|---|---|---|---|---|
| None | 297241 | 0.6359041 | 0.6366664 | 1029.805 |
| V1 | 7328 | 0.4138919 | 0.2969432 | 2130.317 |
| V2 | 5822 | 0.3083133 | 0.3689454 | 4541.016 |
| V3 | 4253 | 0.2031507 | 0.2976722 | 6686.414 |
| hV4 | 2058 | 0.0442177 | 0.0782313 | 7286.004 |
| VO1 | 946 | 0.0073996 | 0.2082452 | 6950.166 |
| VO2 | 830 | 0.0445783 | 0.5084337 | 5926.769 |
| LO1 | 1741 | 0.0063182 | 0.2004595 | 7886.712 |
| LO2 | 1067 | 0.0037488 | 0.3683224 | 6219.818 |
| TO1 | 1354 | 0.0901034 | 0.1964549 | 5826.597 |
| TO2 | 1010 | 0.3356436 | 0.3554455 | 3329.027 |
| V3b | 1006 | 0.1998012 | 0.4423459 | 3825.895 |
| V3a | 3028 | 0.0465654 | 0.3738441 | 5609.348 |
Now get rid of the noisy voxels.
myfile_denoised=myfile[myfile$noisesep=="signal",]
myfile_denoised=na.omit(myfile_denoised)
myfile_denoised=myfile_denoised[myfile_denoised$ROI!="None",]
myfile_denoised=myfile_denoised[myfile_denoised$onscreen==1,]
Show visual field coverage, but this time confined to the screen.
Also plot as a function of ROI.
cov1=ggplot()+geom_hline(yintercept=0,colour='white')+
theme_classic()+theme(panel.background = element_rect(fill='black'))+geom_vline(xintercept=0,colour='white')+geom_circle(alpha=.1,data=myfile_denoised,mapping=aes(x0=X,y0=Y,r=sigma/2,fill=ROI,alpha=.1),colour='transparent')+xlim(-10,10)+ylim(-10,10)+coord_fixed()+facet_wrap(~ROI,nrow=facetval)+ theme(legend.position="bottom")+theme(legend.position = "None")+scale_colour_viridis(discrete = TRUE,option='plasma')
cov1
ggplot(myfile_denoised, aes(x = X, y = Y))+xlim(-10,10)+ylim(-10,10)+geom_hline(yintercept=0,colour='white')+geom_vline(xintercept=0,colour='white')+coord_fixed()+theme_classic()+theme(panel.background = element_rect(fill='black')) +geom_density_2d()+facet_wrap(. ~ROI,nrow=facetval) + stat_density_2d(aes(fill = after_stat(level)), geom = "polygon")+ scale_fill_viridis(option='plasma')
R2plot=ggplot(myfile_denoised,aes(x = reorder(ROI, R2, FUN = mean),y=R2))+theme_classic(base_size=24)+stat_summary(fun = mean, geom = "bar", width = 0.5,aes(fill=..y..))+stat_summary(fun.data = mean_se, geom = "errorbar", width = 0.3)+ theme(legend.position="bottom")+xlab("ROI")+ylab("Variance explained")+ scale_fill_viridis_c(option='plasma',begin=0,end = 0.8)+theme(legend.position = "None")
plot_grid(R2plot,R2c,nrow = 2,rel_heights = c(0.4,0.6))
Splot=ggplot(myfile_denoised,aes(x = reorder(ROI, sigma, FUN = mean),y=sigma))+theme_classic(base_size=24)+stat_summary(fun = mean, geom = "bar", width = 0.5,aes(fill=..y..))+stat_summary(fun.data = mean_se, geom = "errorbar", width = 0.3)+ theme(legend.position="none")+xlab("ROI")+ylab("Sigma")+ scale_fill_viridis_c(option='plasma',begin=0,end=0.5)
plot_grid(Splot,Sc,nrow = 2,rel_heights = c(0.4,0.6))
Eplot=ggplot(myfile_denoised,aes(x = reorder(ROI, ecc, FUN = mean),y=ecc))+theme_classic(base_size=24)+stat_summary(fun = mean, geom = "bar", width = 0.5,aes(fill=..y..))+stat_summary(fun.data = mean_se, geom = "errorbar", width = 0.3)+ theme(legend.position="none")+xlab("ROI")+ylab("Eccentricity")+ scale_fill_viridis_c(option='plasma',begin=0,end=0.3)
plot_grid(Eplot,Ec,nrow = 2,rel_heights = c(0.4,0.6))
SEplot=ggplot(myfile_denoised,aes(x = ecc,y=sigma))+recipe+geom_point(aes(color=ROI))+ theme(legend.position="bottom")+theme_classic(base_size=24)+geom_smooth(method='lm',colour='black')+facet_wrap(~ROI,nrow=facetval)+xlab("Size")+ylab("Eccentricity")+scale_colour_viridis(discrete = TRUE,option='plasma')+ theme(legend.position="none")
SEplot
## `geom_smooth()` using formula 'y ~ x'
polplot=ggplot(myfile_denoised, aes(x = ang)) +
geom_histogram(binwidth = 15, aes(fill=..x..)) +
scale_x_continuous(breaks = seq(0, 360, 60))+theme_classic(base_size=24)+
coord_polar(start=deg2rad(0),direction=1)+facet_wrap(~ROI,nrow=2)+ theme(legend.position="bottom")+ scale_fill_viridis_c(option='plasma')+theme(legend.position = "None")
plot_grid(polplot,Ac,nrow = 2,rel_heights = c(0.5,0.6))